Efficient L1SPIRiT Reconstruction (ESPIRiT) for Highly Accelerated 3D Volumetric MRI with Parallel Imaging and Compressed Sensing

Introduction: Highly accelerated data acquisition is demanded for 3D volumetric MRI. In recent years, many approaches [1,2,3] have been developed to integrate parallel imaging (PI) and compressed sensing (CS) to achieve higher acceleration than either method alone. Among such approaches, L1SPIRiT [3] synergistically combines PI and CS and has proven promising in clinical evaluations. However, this iterative solver is highly computationally intensive and poses difficulty for commonly available platforms. This work was aimed at developing an efficient L1SPIRiT scheme (ESPIRiT) to address this computation challenge. Theory: L1SPIRiT is an iterative algorithm performing PI and CS operations serially in each iteration . The PI operator resynthesizes k-space using a GRAPPA-like convolution kernel (Gk) . This operation can be performed more efficiently with image-domain multiplications : Xn+1(x,y)=Xn (x,y)·GI(x,y) (1), where Xn,n+1(x,y) are temporary image-domain solutions at pixel (x,y), GI is image-domain unliasing coil weights (GI=F(Gk)). The CS operator transforms multi-coil images to sparse domain (w) using wavelet (Ψ) and pursues min ||w||1 using softthresholding (T). The computation of the PI and CS operators is O(NX·NC·Nit) and O(NX·NC·Nit), respectively, where NX, NC·and Nit are the numbers of pixels to reconstruct, coil channels and iterations in the entire reconstruction, respectively. This work intended to reduce computation from the following three perspectives: 1. modified L1SPIRiT to remove NC: The PI operator utilizes k-space correlations and ideally should converge to the “truth” image: X=M·C, where M and C represent spin density and coil sensitivity distributions, respectively. By rewriting equation (1) ((x,y) omitted below for simplicity), we have M·C=M·C·GI (2). By eliminating the common scalar M in (2), we get C=C·GI (3), which means C (size: NC×1) corresponds to the eigenvector of GI (size: NC× NC) with eigenvalue=1 at each pixel. (3) offers an approach to estimate C from GI (Fig. 1A), with which we perform PI & CS in an alternative way. Our PI operator pursues a new solution that is consistent with coil weighting and meanwhile is L2-closest to the previous solution: min ||Xn+1-Xn||2, s.t. Xn+1=M·C. The derived optimal solution is: Xn+1=C·Xn·C / ||C||2 (4). Let Cs=C/||C||2, we can rewrite and split (4) to two operators: S1: Mn+1=Cs·Xn and S2: Xn+1=M n+1·Cs. S1/2 reduces the matrix operation (O(NC)) in (1) to much faster vector operation (O(NC)). Furthermore, an intermediate magnetization image Mn+1 is produced such that CS can now be performed in the coil-combined image pursuing joint sparsity rather than coil by coil. This further reduces the computation of the CS operator by NC×. Additionally, (3) is well-conditioned only at pixels with signals, while in air, (3) produces eigenvalues largely different from 1 (Fig. 1A). Thus, the eigenvalue map of GI can be used to generate an image support (IS) that can eliminate artifacts in air and improve the conditioning of L1SPIRiT . 2. pixel-specific convergence to reduce NX: It is observed that L1SPIRiT convergence is highly pixel-specific (Fig. 2). For most pixels, only a small number of iterations are needed. Taking advantage of this feature, converged pixels can be “checked out” and excluded in later iterations. This can rapidly reduce NX remaining in reconstruction (Fig. 2), which can accelerate S1/2 operators and Fourier and wavelet transforms performed on an increasingly sparser image. 3. PI initialization to reduce Nit: It has been shown that PI can improve the initial condition for L1SPIRiT than conventionally used zero-filling and therefore reduce Nit needed . Accordingly, Poisson-disk k-space sampling (PDS) is replaced by tiled-PDS (tPDS) for efficient PI initialization without sacrificing image quality .